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Introduction 



In the recent years, there has been a growing interest in the study of nonequilibrium sys- 
tems. The states of systems at thermal equilibrium, are a priori given by the Boltzmann- 
Gibbs distribution i.e. there is an energy function E{C) associated to every possible 
configuration C of the system and each configuration C has a weight proportional to 
ex-p[—E{C)/kBT], where ks is Boltzmann's constant and T is the temperature of the 
system. Then all the thermodynamic quantities are obtain by averaging over all the config- 
urations with respective weights. However, in nature there exist a wide variety of systems, 
which are not in thermal equilibrium. The probabilities of different states of these systems 
are not given by the Gibbs distribution, but are determined by the underlying microscopic 
dynamical processes, and are often hard to determine due to lack of a general framework. 

An important class consists of nonequilibrium systems, which when driven by slowly 
varying external forcing, relax through avalanche-like dynamics in response to the exter- 
nal perturbations. Examples include sand or rice piles, forest fires, earthquakes, vortices 
in dirty type II superconductors, solid on solid friction, moving of interfaces in random 
media, disordered ferromagnets and many others (see Jensen 1998, chap. 3). Depending 
on the system, the avalanche is characterized by different physical quantities. For exam- 
ple, in sand piles the system is driven by slowly adding sand grains to the system and the 
avalanche is characterized by the number of sand grains displaced after adding a single 
grain or the lifetime of the avalanche. In earthquakes, it is the energy release and in case of 
ferromagnets it is the size of the domain that flips. The avalanches occurs in various sizes 
in a random sequence, and one is generally interested in the distribution of the avalanche 
sizes. 

What is common in all the systems mentioned above is the existence of threshold and 
multiple metastable states i.e. if the applied external force is less than a critical value the 
system does not response and when the force exceeds the critical value the system passes 
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from one metastable state to another. Due to existence of multiple metastable states, the 
state of a system at a given time depends on the history of evolution (path along which 
the system is evolving in configuration space) and systems exhibit hysteresis in the zero 
frequency limit of external forcing^ 

In this thesis, we study a spin model in the presence of disorder, called random field 
Ising model, introduced by Sethna et al. (1993) in the context of Barkhausen noise and 
hysteresis in disordered ferromagnets. In this model, as the external field increases, the 
magnetization increases as groups of spins flip up together. The dynamics is governed 
by the existence of many metastable states, with large energy barriers separating different 
metastable states. We hope that this study of non-equilibrium response in this model would 
help in the more general problem of understanding the statistical mechanics of metastable 
states in glassy systems. 

The remainder of this chapter is organized as follows, section 1.1 contains a brief re- 
view of theoretical studies of hysteresis in ferromagnets. In section 1.2, we briefly discuss 
Barkhausen effect. In section 1.3, we define the random field Ising model with zero tem- 
perature dynamics, and discuss some earlier results. In section 1.4, we discuss some of 
the equilibrium properties of random field Ising model. Section 1.5 gives an outline of the 
remaining chapters. 



1.1 Hysteresis in ferromagnets 

The studies of hysteresis in magnetic materials has been there in various branches of sci- 
ence, for a long time (see Bertotti 1998). Apart from the intellectual interest, it also has 
wide range of technological applications, from designing transformer cores to memory 
devices. 

Physicists have been looking for a convincing general theory to interpret the phe- 
nomenon of hysteresis in magnets since the time of Rayleigh (1887), who gave the first 
phenomenological theory where the experimental magnetization curves at small field were 
approximated by parabolas. Starting from the demagnetization state (zero magnetization 
in the absence of external field), the magnetizations M± at small fields ±h, are expressed 

t Systems also show hysteresis under periodic forcing. For example, when a ferromagnet is placed in 
oscillation field, the magnetization lags behind its instantaneous equilibrium value and gives rise to hysteresis 
loop. But the area of the hysteresis loop tends to zero in the zero frequency limit of the driving field (Dhar 
and Thomas 1992) 
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Figure 1.1: Rayleigh hysteresis loop 



as 

M± = ah±bh\ (1.1a) 

and when the field is cycled between small zti^max, the lower and the upper curves of the 
hysteresis loop are represented by 

M± = (a + bH,^,,)h ± h{h' - HlJ. (1.1b) 

In Fig. 1.1, we have shown the magnetization curves starting from the demagnetization 
state given by Eq. (1.1a) and the Rayleigh hysteresis loop given by Eq. (1.1b). 

In Weiss's (1907) theory of ferromagnetism, he postulated the existence of a powerful 
internal "molecular field" in ferromagnet materials, which would tend to tries to align the 
magnetic moments along one direction. It agrees with some of the experimental cases 
where it is possible to attain a large saturation magnetization by the application of a very 
weak magnetic field [see Fig. 1.2]. However, it did not explain the fact that, it is also 
possible for the magnetization to be zero (or nearly zero) in the absence of a magnetic field. 
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Figure 1.2: Magnetization curve of single crystal of silicon iron. (Williams and Shockley 1949). 



Weiss further made an assumption that, a ferromagnet can be subdivided into regions called 
magnetic domains. In each individual domain, the the magnetic moments are aligned along 
the molecular field, but the orientation of the spontaneous magnetizations in each domain 
distributed randomly inside the sample and hence the resultant magnetization could be zero 
in the absence of a external field, even at very low temperature. If a external field is applied 
opposite to the magnetization direction, a domain reverses the direction of magnetization 
when the external field exceeds a critical value He. Therefore, if a gradually increasing 
external field is applied, domains whose magnetization vectors are at an angle (tt — 9) with 
the external field, will suddenly reverse direction when the external field exceeds He / cos 9. 
This results a finite bulk magnetization for external field h > He. A comprehensive review 
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of the physical principles of domain theory and of the crucial experiments which bear 
directly on the foundation of the subject, may be found in an article by Kittel (1949). 

Preisach (1935) introduced a modified domain model, in which he assumed that mate- 
rial is composed of many small domains and each of them possesses a rectangular hystere- 
sis loop. The interaction between domains are represented by a local field acting on each 
domain. Thus each domain has two different coercive fields a and (3 for the increasing 
and decreasing branches respectively. The ensemble of domains is then described by the 
distribution function P(a, /3) of the values of a and (3 and hysteresis loops are obtained by 
taking the weighted sum of magnetization in all the domains. 

Sethna et al. (1993) proposed the random field Ising model with the zero temperature 
dynamics as a simple theoretical model for the Barkhausen noise and hysteresis in dis- 
ordered ferromagnets. In this model, magnetic domains are represented by Ising spins 
{si = ±1) and the external field is coupled to these spin. In contrast to the Preisach 
model of hysteresis, where interactions between the individual hysteresis units (domains) 
are ignored, in the random field Ising model the spins interact ferromagnetically with their 
neighbors. The homogeneities and disorder in materials are modeled by introducing a un- 
correlated random field acting on each domain, chosen at random from some distribution. 
Since the domains interact ferromagnetically, flipping of a domain at some external field 
may force the neighboring domains to flip as well in the same direction, thus leading to an 
avalanche of domain flips, which is analogue of a Barkhausen pulse in real magnets (for a 
comparison, see Fig. 1.5 and Fig. 1.6). 

In the models discussed above, the hysteresis does not depends on the rate at which the 
external field is varied i.e. relaxation time from one metastable state to another is much 
larger to the rate at which the system is driving. In contrast, there are also other models 
studied in the context of rate-dependent hysteresis, where the system exhibits hysteresis 
only when it is driven at a finite rate (see Chakrabarti and Acharyya 1999, for a recent 
review). Hysteresis in the A^-vector model was widely studied by many authors (Rao 
et al. 1990, Dhar and Thomas 1992, 1993, Somoza and Desai 1993). It was shown that 
in all dimensions d > 2, for > 2 at low frequency cu and low amplitude Hq of the 
driving field the area of the hysteresis loop scales as [Houyl'^ with logarithmic corrections. 
At high frequencies the area varies as Hq/lo. For any Hq, there is a dynamical phase 
transition separating these two frequency regimes. Above the critical frequency u{Ho), the 
hysteresis loop does not posses inversion symmetry. Using the nucleation theory Dhar and 
Thomas (1993) showed that for = 1 and o? > 1, the area of the hysteresis loop scales as 
|rin(i/ocu)|-i/('^-i) forcu < Hq. 
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1.2 Barkhausen effect 

The first evidence for the existence of ferromagnetic domains was from an experiment 
by Barkhausen (1919). His experiment consists in amplifying the voltage induced in a 
secondary pick-up coil wound around a ferromagnetic sample, while the sample is being 
magnetized by a continuous variation of external magnetic field [Fig. 1.3]. He observed a 
noise induced in the pick-up coil, corresponds to a sudden, discontinuous jumps in magne- 
tization [Fig. 1.4]. These jumps are interpreted as discrete changes in the size or rotation 
of ferromagnetic domains. An elementary introduction of the Barkhausen effect may be 
found in the textbook by Feynman et al. (1977). 




Figure 1.3: Barkhausen effect 



In the recent years, there has been a great interest in the study of the statistical prop- 
erties of the Barkhausen noise. A typical train of barkhausen noise signals observed in 
experiments is shown in Fig. 1.5. Three basic physical quantities that describe a single 
Barkhausen noise signal in an experiment, are signal duration, area of the signal and the 
energy released during the signal occurrence. It is observed that distribution of these quan- 
tities follow power law over a few decades with a cut-off as shown in Fig. 1.7 (see Spaso- 
jevic et al. 1996, and references therein). This power law tail in the Barkhausen avalanche 
distribution was interpreted by Cote and Meisel (1991) as an example of self-organized 



Chapter 1. Introduction 



7 




Figure 1.4: Barkhausen jumps 



criticality^ . But Perkovic et al. (1995) have argued that large bursts are exponentially rare, 
and the approximate power-law tail of the observed distribution comes from crossover ef- 
fects due to nearness of a critical point. 

Barkhausen effect is also widely used as a noninvasive material characterization tech- 
nique for ferromagnetic materials (see Sipahi 1994, for an overview). 



1.3 Hysteresis in random field Ising model 

The nonequilibrium random field Ising model was proposed by Sethna et al. (1993) as a 
model for Barkhausen noise and hysteresis in ferromagnets. The model is defined on a 
lattice. At each lattice site i, there is a Ising spin Si = ±1, which interacts with nearest 
neighbors through a ferromagnetic exchange interaction (J > 0). Spins {si} are coupled 
to the on-site quenched random magnetic field hi and the external field h. The Hamiltonian 
of the system is given by 

i7 = -J^ SjSj - ^/ijS— /i^Sj, (1.2) 

where (z, j) denotes that the sum runs over nearest neighbor pairs of spins on sites i and 
j. We assume that {hi} are quenched independent identically distributed random variables 

^In self-organized critically, systems exhibit critical behavior (power law correlations), without fine tun- 
ing any parameter (for an overview, see: Dhar 1999, Jensen 1998, Bak 1997). 
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Figure 1.5: An example of experimental Barkhausen signal (voltage pulse produce from a pickup 
coil around a fen^omagnet subjected to a slowly varying applied field) (Urbach et al. 1995). 




Figure 1.6: Time series of the avalanches (the number of spin flips at a given field) in the random 
field Ising model on a square lattice of size 200 X 200. From one avalanche to the next avalanche is 
considered as one time step. 
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Figure 1.7: Experimental data for the distribution of Barkhausen signal durations, areas and ener- 
gies (Spasojevic et al. 1996). 
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with the probability that the value of the random field at site i lies between hi and hi + dhi 
being (f>{hi) dhi. 

When the external field is changed, the system relaxes to a stable spin configuration 
through zero-temperature Glauber single- spin- flip dynamics (Kawasaki 1972, see), which 
is specified by the transition rates 



where AE is the change of energy in the system as a result of the spin flip. Therefore, a 
spin-flip is allowed only if the process lowers energy. We assume that the external field 
is increased adiabatically, i.e. F ^ cu, the rate at which the magnetic field h is increased. 
Thus if the spin-flip is allowed, it relax instantly, so that the spin Si in a stable configuration 
is parallel to the net local field at the site: 



Note that the limit uj/T — > is taken after the limit T — > 0. If the limits are taken in 
the reverse order, the state of the the system at each h, is the equilibrium state for all finite 
T and the hysteresis loop area goes to zero. 

We start with h = —oo, when all spins are down and slowly increase h. As we increase 
h, some sites where the quenched random field is large positive will find the net local 
field positive, and the spin at that site will flip up. Flipping a spin makes the local field 
at neighboring sites increase, and in turn may cause them to flip. Thus, the spins flip in 
clusters of variable sizes. If increasing h hy a very small amount causes s spins to flip 
up together, we shall call this event an avalanche of size s. As the applied field increases, 
more and more spins flip up until eventually all spins are up, and further increase in h has 
no effect. 

As an illustrative example, consider a four by four square lattice with periodic bound- 
ary condition and a particular realization of quenched random fields, which is shown in 
Fig. 1.8(a). We set J = 1. Now start with h = — oo, when all spins are down [Fig. 1.8(b)] 
and slowly increase it. A spin with m up neighbors, flips up at h, if the quenched random 
field at the particular site hi > i — 2m — h. Therefore, when the external field just exceeds 
the value 1.1, the local field at the site where hi = 2.9, becomes positive and the spin at 
that site flips up. This increases the local fields at its neighboring sites by 2 J and as a result 
spins at some of these sites flip up [Fig. 1.8(c)]. These process continues till there is no 
more sites where the local field is positive at that external field. In the figure, we denote 




(1.3) 
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Figure 1.8: (a) Quenched random fields at various sites, (b) - (f) Stable spin configurations at 
various external fields. The black spins ai^e inactive spins at that particular field and the colored 
spins are part of the avalanche. The colors specify the order in which spins flip during the avalanche, 
(g) Clusters of spins, which flip during one avalanche, (h) Color map showing the order of events, 
(i) Magnetization curve, corresponding to evolution. 
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active spins with different color according to order at which they flip during the avalanche. 
As shown in Fig. 1.8(b) - (f), the system passes from one stable configuration to another, 
as the external field is increased, till all the spins in the system flip up. In Fig. 1.8(g), we 
show the different clusters of spins, which flip at different fields. Fig. 1.8(i) shows the 
corresponding magnetization. 

Sethna et al. (1993) studied the model with the infinite-range interaction (mean field 
theory), where every spin is coupled to all N other spins with coupling J/N. They found 
that there exists a critical value of disorder Ac (which in the case of a Gaussian distribution 
of random fields is = ^J{2/tt)J ), below which the hysteresis curve displays a jump due 
to an infinite avalanche of spin flips, which spans the system. Above the critical disorder 
systems show smooth magnetization curve without macroscopic jumps. However, this 
mean field theory does not show any hysteresis for disorder A > Ac. Dahmen and Sethna 
(1993, 1996) studied the hysteresis loop critical exponents expanding about mean field 
theory in 6 — e dimensions. A power- low distribution with avalanche of all sizes is seen only 
at the critical value of the disorder. However, the numerical simulations by Perkovic et al. 
(1995) indicate that the critical region is remarkably large: almost three decades of power- 
law scaling in the avalanche size distribution remain when measured 40% away from the 
critical point. Therefore, they argued that several decades of scaling seen in experiments 
need not be self-organized criticality, as many of the samples might have disorders within 
40% of the critical value. 

Interestingly the model can be solved exactly on a Bethe lattice for the magnetization 
on the hysteresis curve for arbitrary distribution of random fields (Dhar et al. 1997). In 
contrast to the infinite-ranged mean field theory, the calculation on Bethe lattice shows 
hysteresis even for large disorder. Another interesting result of the Bethe lattice calculation 
is that, the first order jump in the magnetization disappear for coordination number of the 
lattice less than 4. Only for coordination number 4 and above, there exists a critical value 
of disorder below which there is a jump discontinuity in the magnetization. 



1.4 Equilibrium properties of random field Ising model 

In this thesis we are interested in the nonequilibrium properties of the random field Ising 
model. However, it is useful to recall the equilibrium properties of this model, which 
has been an important problem in statistical physics for a long time. This model has a 
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number of interesting realization in nature. A recent review of earlier work on this model 
may be found in an article by Nattermann (1998). This model was first studied by Imry 
and Ma (1975), in the context of possible destruction of long-range order by arbitrarily 
weak quenched disorder. The pure Ising model with nearest neighbor interaction Hq = 
— J J2(i,j) ■Si'Sj , is known to have a ferromagnetically ordered phase in all dimensions d > 1. 
When the random field term — ^ hiSi is introduced, it acts against the order. Imry and Ma 
(1975) argued that arbitrarily weak disorder destroys long-ranged ferromagnetic order in 
dimensions d <2. The argument goes as follows: 




Figure 1.9: Domain of reverse spins. 



If we consider a domain of reverse spins [Fig. 1 .9], of linear size ~ L, the domain wall 
energy is ~ L'^'^ . However, according to the central limit theorem, if the random field 
has short-range spatial correlations, the fluctuation in the magnetic field energy in such 
domains is ~ L''/^. Thus, by splitting into domains of size L, the system will gain bulk 
energy of 0{L'^l'^) per domain, and loose a surface energy a surface energy of 0{L'^~^) 
per domain. Thus, whenever d < 2, there will exists a large enough L, for an arbitrarily 
small random field, where it will become energetically favorable to the system to break 
into domains of that size. 

The argument by Imry and Ma (1975) suggests that the lower critical dimension^ is 
di > 2, rather than di = 2, because other mechanisms could destroy long-range order 
in higher dimensions. It is widely believed that the upper critical dimension^ is c/„ = 6, 
instead of c/„ = 4 for the pure Ising system. However, whether the lower critical dimension 
di = 2 or 0?; = 3, was a matter of a long controversy, but has now been established that 

^The dimension below which long-range ferromagnetic order cannot exist. 

+The dimension above which the critical exponents are those of the Gaussian fixed point. 
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di = 2. Imbrie (1984) showed that if the disorder is small, the model in dimension o? = 3 
exhibits long-range order at zero temperature. Aizenmann and Wehr (1989) rigorously 
proved uniqueness of the Gibbs state in o? = 2, i.e. absence of any phase transition, in 
agreement with the Imry-Ma prediction. 

As far as an exact calculation of thermodynamic quantities is concerned, there are only 
a few results. For example, Bruinsma (1984) studied the random field Ising model on a 
Bethe lattice in the absence of an external field and for a bivariate random field distribu- 
tion. There are no known exact results for the average free energy or magnetization, for a 
continuous distribution of random field, even at zero temperature and in zero applied field. 



1 .5 Outline of this tliesis 

In this thesis, we study the nonequilibrium ferromagnetic random field Ising model with 
zero temperature Glauber single flip dynamics. The remaining chapters of the thesis are 
organized as follows: 

In chapter 2, we discuss some special properties of the model that makes the analytical 
treatments possible. We briefly recapitulate the derivation of self-consistent equations for 
the magnetization in the model. 

In chapter 3, we use a similar argument to construct the generating function for the 
avalanche distribution for arbitrary distribution of the quenched random field. In sec- 
tion 3.2, we consider the special case of a rectangular distribution of the random field. 
In this case, we explicitly calculate the probability distribution of avalanches, for the for 
Bethe lattices with coordination numbers z = 2 and 3. In section 3.3, we analyse the self- 
consistent equations to determine the form of the avalanche distribution for some general 
unimodal continuous distributions of the random field. In chapter 4, we derive the self- 
consistent equations for the magnetization on minor hysteresis loops on a Bethe lattice, 
when the external field is varying cyclically with decreasing magnitudes. We also discuss 
some properties of stable configurations, when the external field is varying. 

In chapter 5, we study the model with an asymmetric distribution of quenched fields, 
in the limit of low disorder in two and three dimensions. We relate the spin flip process to 
bootstrap percolation, and find nontrivial dependence of the coercive field on the coordina- 
tion number of the lattice. 

Chapter 6 contains a discussion of our results, and some concluding remarks. 
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Some algebraic details of the analytical solution for the distribution of avalanche sizes, 
for the rectangular distribution of quenched fields on Bethe lattice are relegated to ap- 
pendix A. 

Part of this work has appeared in journals as refereed papers. Though it is mostly a 
repetition of the material presented in chapter 3 and chapter 5, for the convenience of the 
reader, we have reproduced these papers as an appendix (reprints) to this thesis. 



Chapter 2 



Earlier exact results on hysteresis 
in random field Ising model 



The difficulty of solving various mathematical equations describing actual physical situa- 
tions leads to various approximation methods. These approximation method can be clas- 
sified into two categories: One in which the approximation is made in the mathematical 
equations itself and another in which the physical system is simplified. Into the second 
category fall many lattice model systems. Again in higher dimensions the lattices contain 
closed circuits which makes the model difficult to solve. Thus, one considers the problem 
in the mean-field theory, where the underlying lattice structure becomes irrelevant or on a 
different lattice where it can be solved exactly. Bethe lattice or Cayley tree is one in which 
there is no circuits at all which makes the model easier to solve. The simplicity of the 
lattice motivates one to study various systems on a Bethe lattice. 

In this chapter we briefly discuss derivation of hysteresis curve in the random field 
Ising model in the mean field theory (infinite-range interaction) and on the Bethe lattice. 
In sections 2.1 and 2.2 we discuss two properties of zero temperature random field Ising 
model, namely the return point memory effect and the abelian nature of spin-flips, which is 
used to set up self-consistent equation to determine magnetization on the Bethe lattice and 
later in other chapters. 



2.1 Return point memory 

Sethna et al. (1993) showed that the RFIM exhibits the following return point memory 
effect: Suppose we start with h = — oo, and all spins down at t = 0. Now we change 
the field slowly with time, in such a way that h{t) < h{T), for all times t < T. Then 
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(a) 



(b) 



Figure 2.1: (a) change of external field with time, (b) magnetization curve vs. external field given 
by (a). When the external field returns to the previous extremal value, the magnetization returns to 
the value at that field i.e. M{h{ti)) = M{h{t2)). 

the configuration of spins at the final instant t = T does not depend on the detailed time 
dependence of h{t), and is the same for all histories, so long as the condition h{t) < h{T) 
for all earlier times is obeyed. In particular, if the maximum value h{T) of the field was 
reached at an earlier time ti, then the configuration (and hence the magnetization) at time 
T is exactly the same as that at time ti [Fig. 2.1]. 

Consider two spin configurations C{si, S2, • • • , Sn} and C"{s[^ s^^ . . . , s'^}. If Si > s'- 
for each site i, the configurations C and C" are called partially ordered, C > C". Let two 
configurations C{t) and C"{t) be evolve under the field h{t) and h'{t) respectively. Suppose 
the initial configurations C(0) and C"(0) are partially ordered such that C(0) > C"(0) 
and the fields satisfy h{t) > h'{t). Then if a spin is up in configuration C'{t), the 
corresponding spin Si{t) in configuration C{t) must be up, since the local field li{t) in 
C{t) can not be less than in C'{t). Therefore the configurations C{t) and C'{t) will 
always remain partially ordered, C{t) > C"{t). This is the no passing property of the 
system. An earlier treatment of "no passing" rule was given by Middleton (1992) in the 
context of charged-density waves. 

Let us consider the Fig. 2.2. The configuration A is reached by increasing the field from 
a lower value to hi. On increasing the field monotonically from hi to h2, configuration B 
is reached. Naturally, the configurations A and B are partially ordered such that 



B> A. 



(2.1a) 
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Similarly when the field is decreased monotonically from h2 to hi, partial ordering exists 
between the configuration B and the final configuration C such that 

C <B. (2.1b) 

Since during the evolution from A to C, the field h{t) satisfies h{t) > hi, 

C >A. (2.1c) 

Now suppose the configuration C evolve to a configuration D when the field is again 
increased monotonically from hitoh2. Since the partial ordering is preserved by dynamics, 
fromEq. (2.1b), 

D <B, (2.1d) 

and from Eq. (2.1c), 

D>B, (2.1e) 

as A evolves to B, when the field is increased from hi to h2 . From Eq. (2. Id) and Eq. (2. le), 
we must have 

D = B, (2.1f) 
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i.e. the system returns exactly to the same earlier configuration, when the field is decreased 
from a maximum value h2 and then increased to the same value. The same memory effect 
extends to subcycles within the cycles and so on. 



2.2 Abelian property 

Because of the previous property, we may choose to increase the field suddenly from — oo 
to /i(r) in a single step. Then, once the field becomes h = h{T), several spins would have 
positive local fields. Suppose there are two or more such flippable sites. Then flipping any 
one of them up can only increase the local field at other unstable sites, as all couplings are 
ferromagnetic. Thus to reach a stable configuration, all such spins have to be flipped, and 
the final stable configuration reached is the same, and independent of the order in which 
various spins are relaxed. This is the abelian property of relaxation (Dhar et al. 1997). 
Using the symmetry between up and down spins, it is easy to see that the abelian property 
also holds whether the new value of field h" is greater or less than its initial value h' so 
long as one considers transition from a stable configuration at h' to a stable configuration 
at h". 



2.3 Hysteresis in tlie infinite-range interaction model 

In this section, we will briefly discuss the results obtained by Sethna et al. (1993), on the 
hysteresis in the random field Ising model with infinite-range interaction. In this mean field 
theory, every spin is coupled to all N other spins with coupling J jN. The Hamiltonian is 
given by 

H = -j^w2si\ -h^s,-^h,s,. (2.2a) 

\ i / i i 

Now the interaction of a spin with other spins is replaced by its interaction with the mag- 
netization M{h) of the system. The Hamiltonian then takes the form 

H = - J2{JM + h + hi)s,, (2.2b) 

i 

i.e. the effective local field at site is JM + h + hi. The spin at this site will flip up if 
this field is positive i.e. the quenched random field hi at this site exceeds —.JM — h. This 
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happens with probability 

/oo 
(l){h,) dh,. 
-JM(h)-h 

Therefore, the average magnetization satisfies the self-consistent equation 

/oo 
^{hi)dh,-l. (2.3) 
-JM{h)-h 

Note that, for symmetric distributions of random fields, M(0) = is the trivial solution at 
/i = 0. Now, if M(0) = is the only solution at /i = 0, then there is no hysteresis. To 
have other nontrivial solutions for M(0), the slope of the expression on the right hand side 
of Eq. (2.3), as a function of M(0) must be greater than unity at M(0) = 0. At /i = 
and near M(0) = 0, the right hand side of Eq. (2.3) can be approximated as 2M(0) J</)(0). 
Therefore, the condition that the Eq. (2.3) to has multi-valued solution is 

m > ^j- (2.4) 

This condition corresponds to a critical disorder strength Ac (width of the random field 
distribution), above which there is no hysteresis i.e. the magnetization follows the same 
curve in the increasing and decreasing field. Below Ac, the magnetization curves in the 
increasing and decreasing field are different near h = 0, i.e. the system exhibits hysteresis. 
Moreover, there is a critical field, where the magnetization jumps from one solution to 
another one. In a specific case, where the random field distribution is Gaussian, 

using the condition given by Eq. (2.4), the critical value of disorder is obtained as 

Ac = \/{2/tt)J. (2.6) 

Figure 2.3 shows the magnetization curves for this mean-field at various values of disorder 
A < Ac, A = Ac and A > Ac for Gaussian distribution of random fields. Note that 
hysteresis and jump in the magnetization exist only below a critical disorder [Fig. 2.5(c)]. 

Sethna et al. have studied in detail the case of critical disorder, and the power-law 
divergence of various quantities at this critical point. The special value of disorder does 
not seem to be particularly important and we shall not discuss it here. 



2.4 Hysteresis on the Bethe lattice 

The shortcoming of the treatment discussed in the previous section is that the pair couplings 
are weak and no correlations and short-range order. One can keep mean field theory, but 
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h/J 



(c) 

Figure 2.3: Magnetization curve for the random field Ising model with infinite-range interaction 
at various values of disorder: (a) A = J, (b) A = Ac = ^/ (2/7r)J and (c)A = 0.5J, for the 
Gaussian random field distribution given by Eq. (2.5). The dashed line in (c) shows the third root 
of the self-consistent equation for magnetization, given by Eq. (2.3). 
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add correlations by working on a Bethe Lattice. 

The advantage of working on the Bethe lattice is that the usual BBGKY hierarchy of 
equations for correlation functions closes, and one can hope to set up exact self-consistent 
equations for the correlation functions. The fact that Bethe's self-consistent approximation 
becomes exact on the Bethe lattice is useful as it ensures that the approximation will not vi- 
olate any general theorems, e.g. the convexity of thermodynamic functions, sum rules etc. 
In the presence of disorder, in spite of the closure of the BBGKY hierarchy, the Bethe ap- 
proximation is still very difficult, as the self-consistent equations become functional equa- 
tions for the probability distribution of the effective field. These are not easy to solve, and 
available analytical results in this direction are mostly restricted to one dimension, or to 
models with infinite-ranged interactions. On the Bethe lattice, for short-ranged interac- 
tions with quenched disorder, e.g. in the prototypical case of the ±J random-exchange 
Ising model, the average free energy is trivially determined in the high temperature phase, 
but not in the low-temperature phase. It has not been possible so far to determine even the 
ground- state energy exactly despite several attempts. 



Figure 2.4: A Cayley tree of coordination number 3 and 4 generations. 



The random field Ising model model on a Bethe lattice is special in that the zero- 
temperature nonequilibrium response to a slowly varying magnetic field can be determined 
exactly (Dhar et al. 1997). To be precise, the average non-equilibrium magnetization in 
this model can be determined if the magnetic field is increased very slowly, from — oo to 
+00, in the limit of zero temperature. It thus provides a good pedagogical model to study 
the slow relaxation to equilibrium in glassy systems. 
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The usual way to solve a problem on Bethe lattice is to consider the problem on a 
Cayley tree and calculate all the thermodynamic quantities in the deep inside of the tree. 
Consider a uniform Cayley tree of n generations where each non-boundary site has a coor- 
dination number z and the boundary sites have coordination number 1 [see Fig. 2.4]. The 
first generation consists of a single vertex. The r-th generation has z{z — 1)''"^ vertices for 
r > 2. At each vertex there is a Ising spin. 

Because of the return point memory, to find the magnetization at field h in the lower 
half of the hysteresis loop, we start with h = —oo, when all spins are down and increase 
the field to /i in a single step. Now at this field, since the spins can be relaxed in any order 
(abelian property), we relax them in this: First all the spins at generation n (the leaf nodes) 
are relaxed. Then spins at generation n — 1 are examined, and if any has a positive local 
field, it is flipped. Then we examine the spins at generation n — 2, and so on. If any spin 
is flipped, its descendent are reexamined for possible flips^ In this process, clearly the 
fiippings of different spins of the same generation r are independent events. 

Let P^'''\h) be the probability that a spin on the r-th generation will be flipped when 
its parent spin at generation r — 1 is kept down, the external field is h, and each of its 
descendent spins has been relaxed. As each of the z — 1 direct descendents of a spin is 
independently up with probability p(''+i), the probability that exactly m of them are up is 
— p(''+i)]^-i-'". Suppose we pick a site at random in the tree away from 
the boundary, the probability that the local field at this site is positive, given that exactly m 
of its neighbors are up, is precisely the probability that the local field hi at this site exceeds 
[{z — 2m) J — h]. We denote this probability by Pm{h). Clearly, 

Pm{h) = / (pih,) dh,. (2.7) 

J (z-2m)J-h 

Now it is straightforward to write down a recursion relation for P^*") in terms of p(''+i): 

P^'\h) = E " ^] [P^^^'\h)Y [l - P^^^^\h)\-'-^ p.^{h). (2.8) 

m=0 \ ^ / 

Given a value of h, we can determine pm(^) using Eq. (2.7). Then using Eq. (2.8), and the 
initial condition P(") = po{h), P^''^ can be determined for all r < n. For r ^ n, these 

^This step is not really necessary if we are only interested in determining the magnetization at the site O. 
Skipping this step leads to considerable simplification of the relaxation process: First the spins of generation 
n are examined, then those of (n — 1) etc. till we finally examine the spin at O. No spin is checked more 
than once. The resulting configuration is not fully relaxed, but it is easy to prove that further relaxation will 
not change the state of the spin at O. The argument can be extended to show that the probability that an 
avalanches starting at O is of size s also is the same in this partially relaxed state as in the fuUy relaxed state. 
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probabilities tend to limiting value, lim„^oo P^"^^ = P*, which satisfies the equation 

P'-ih) = E " ^] [P'-ihT [1 - P''{h)r'-'- Pm{h). (2.9) 

m=0 \ / 

This is a polynomial equation in P*{h), which can be solved in terms of {pm{h)}. 
Finally, for the spin at O, there are z downward neighbors, and the probability that it is up 
is given by 

Prob(3o = +l\h)=j2('] [P*ihT [1 - P*ih)Y-"' Pm{h). (2.10) 

m=0 \ / 

Because all spins deep inside the tree are equivalent, Prob(so = +1 1 h) determines the 
average magnetization for all sites deep inside the tree. This determines the lower half of 
the hysteresis loop. The upper half is obtained similarly. 

For the three coordinated {z = 3) Bethe lattice, the self-consistent equation satisfied by 
P* [Eq. (2.9)] is quadratic and from physical arguments at least one root must vary between 
and 1 continuously with h for any value of disorder strength A. Hence the magnetization 
given by Eq. (2.10) is also a continuous function of h. This is also the case with a linear 
chain (z = 2), where the self-consistent equation [Eq. (2.9)] is linear. 

On the other hand, the situation is quite different for z > i. For example, for z = 
Eq. (2.9) is cubic, which has either one or three real roots which will vary with h. Figure 2.5 
shows this variation for two values of disorder of the random field distribution given by 

4>{hi) = ^sech\h,/A). (2.11) 

Note that for large disorder, there is only one real root which vary continuously from 
to 1, giving rise to a continuous magnetization curve as shown in Fig. 2.6. But for small 
disorder, P*{h) as a function of h shows a "S" shaped curve, where at some value of h, 
two real root merge to becomes imaginary and disappear from the real plane. Therefore, 
as we vary h, on the physical ground initially P*{h) takes the lower value till the point 
where it becomes complex and at that point it jumps to the upper value, giving rise to a 
jump discontinuity in the corresponding magnetization curve [Fig. 2.6]. 

This can be generalized to higher coordination number, where the mechanism of two 
real solutions of the polynomial equation Eq. (2.9) merging and both becoming complex is 
still the same. 
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Figure 2.5: Variation of P*{h) with h for the Bethe lattice with z = A, and the random field 
distribution given by Eq. (2. 11). 




Figure 2.6: Magnetization as a function of increasing field for the Bethe lattice with z = A and the 
random field distribution given by Eq. (2. 11). 



Chapter 3 



Distribution of avalanche sizes 
on the Bethe lattice 



In this chapter, we study the distribution of avalanche sizes in the random field Ising model 
on a Bethe lattice (Sabhapandit et al. 2000). This chapter is organized as follows. In 
section 3.1, we set up a self-consistent equation for the generating function Q{x) of the 
probability Qn, that an avalanche propagating in subtree flips exactly n more spins in the 
subtree before stopping, for arbitrary distribution of quenched random fields. Then we 
expressed the generating function G{x) of distribution of avalanche sizes, in terms of Q{x). 
In section 3.2, we consider the special case of a rectangular distribution of the random 
field. In this case, we explicitly solve the self-consistent equations for Bethe lattices with 
coordination numbers z = 2 and 3. However, this case is non-generic. For small strength of 
disorder, the magnetization jumps from —1 to +1 at some value of the field, but for larger 
disorder, when the system shows finite avalanches, there is no jump in magnetization and 
the distribution function decays exponentially for large s. In section 3.3, we analyse the 
self-consistent equations to determine the form of the avalanche distribution for unimodal 
continuous distributions of the random field. We find that for coordination number z > i, 
the magnetization shows a first order jump discontinuity as a function of the applied field at 
some field-strength h^isc, for weak disorder. Just below h = hdisc, the avalanche distribution 
has a universal (—3/2) power-law tail. 



3.1 Generating function for avalanche distribution 

Consider a Cay ley tree rooted at O, of generation [Fig. 2.4]. We will be interested in 
the portion of the tree where generation r ^ A^, in the limit — > oo. Now consider the 
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state of the system at external field h, and all the flippable sites have been flipped. We 
increase the field by a small amount dh till one more site becomes unstable. We would like 
to calculate the probability that this would cause an 'avalanche' of n spin flips. Since all 
sites deep inside are equivalent, we may assume the new susceptible site is the site O. 

It is easy to see that this avalanche propagation is somewhat like propagation of infec- 
tion in the contact process on the Bethe lattice. The 'infection' travels downwards from the 
site O which acts as the initiator of infection. If any site is infected, then it can cause infec- 
tion of some of its descendents. If the descendent spin is already up, it cannot be flipped; 
such sites act as immune sites for the infection process. If the descendent spin is down, it 
can catch infection with a finite probability. Furthermore, this probability does not depend 
on whether the other 'sibling' sites catch infection. Infection of two or more descendents 
of an infected site are uncorrelated events. Thus, we can expect to find the distribution 
of avalanches on the Bethe lattice, as for the size distribution of percolation clusters on a 
Bethe lattice (Stauffer and Aharony 1992). However, a precise description in terms of the 
contact process is complicated, as here the infection spreads in a correlated background of 
'immune' (already up) spins, and the probability that a site catches infection does depend 
on the number of its neighbors that are already up. 

We start with the initial configuration of all spins down. Now increase the external field 
to the value h. Consider a site X at some generation r > 1 of the Cay ley tree [Fig. 3.1]. 
We call the subtree formed by X and its descendents Tx, the subtree rooted at X. We keep 
its parent spin Y at generation r — 1 down, and relax all the sites in Tx at the uniform field 
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h. If X is far away from the boundary, the probability that spin at X is up is P*{h), which 
is obtained by solving the self-consistent equation given by Eq. (2.9). The conditional 
probability that spin at a descendant of X is up, given that the spin at X is down is also 
P*{h). We measure the response of Tx to external perturbation by forcibly flipping the 
spin at Y (whatever the local field there) and see how many spins in this subtree flip in 
response to this perturbation. Let Qn be the probability that the spin at X was down when 
Y was down and n spins on the subtree Tx flip up if Sy is flipped up. Here allowed values 
of n are 0, 1, 2, . . .. Clearly, we have 



^* + E^- = l- (3-1) 

n=0 

We define Q{x) be the generating function of Qn as, 

oo 

Q{x) = E Qnxr (3.2) 

n=0 

Clearly, 

Q{x = 0) = go, (3.3) 

Q{x = l) = l-P\ (3.4) 

It is straight forward to write the self-consistent equation for Q{x). Let us first relax all 
spins on Tx keeping X and Y down. The probability that exactly m the descendents of X 
are turned up in this process be denoted by Pr(m). Clearly 

Pr(m) = r ~ M P*"!! - py-^-^. (3.5) 
\ m J 

For a given m, the conditional probability that local field at X is such that spin remains 
down, even if Y is turned up is 1 — Pm+i ■ Summing over m, and using the expression for 
Pr(m) above, we get 

m=0 \ / 

We can write down an expression for Qi similarly. In this case, if m of the direct 
descendents of X are up when Y is down, the local field at all the remaining z — l — m direct 
descendents must be such that they remain down even if X is flipped up. This probability 

is (^;;^)P*'"go"^"'"- The local quenched field at X must satisfy {z - 2m) J - h > hx > 
[z — 2m — 2) J — h. The probability for this to occur is Pm+i — Pm - Hence we get 

m=0 \ / 
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The equation determines Qn for higher n can be written down similarly. It only involves 
the probabilities Qm with m < n for the descendent spins. Formally we can write 



n \ m I 

m=0 \ / 



oo z — l — m 

n Qn,S{Y^n„n - 1) 

{n,}=0 i=l 



{Pm + 1 ~ Pm )j 



(3.8) 



where 



( 1 for i = j 
[ for i ^ j 
is the Kroneker delta. 

These recursion equations are expressed more simply in terms of the generating func- 
tion Q{x). Multiplying both sides by and then summing over n, we see that the self- 
consistent equation for Q{x) is 



z-l 



m=0 



Q{x) = Q{x = 0) + xJ2[' ^ ) P*" Q{xy-'-^{p^+, - p„ 



(3.9) 



This is a polynomial equation in Q{x) of degree z — 1, whose coefficients are functions 
of h through P*{h) and Pm{h). Using Eq. (2.9) and Eq. (3.6), it is easily checked that for 
X = 1, the ansatz Q{x = 1) = 1 — P* satisfies the equation, as it should. To determine 
Q{x) for any given external field h, we have to first solve the self-consistent equation for 
P* [Eq. (2.9)]. This then determines Q{x = 0) using Eq. (3.6), and then, given P* and 
(5(0), we solve for (5( a;) by solving the (2; — l)-th degree polynomial equation Eq. (3.9). 

Finally, we express the relative frequency of avalanches of various sizes when the ex- 
ternal field is increased from h to h + dh in terms of Q{x). Let Gs{h) dh be the probability 
that avalanche of size s is initiated at O. We also define the generating function G{x\h) as 



G{x\h) = Y,Gs{h)x\ 



(3.10) 



Consider first the calculation of Gs{h) for s = 1. Let the number of descendents of O 
that are up at field hhe m. For the spin at site O to be down at h, but flip up at h + dh, 
the local field ho must satisfy [{z — 2m) J — [h + dh)] < ho < [{z — 2m) J — h]. This 
occurs with probability (f>{zJ — 2m J — h) dh. Each of the [z — m) down neighbors of O 
must not flip up, even when so flips up. The conditional probability of this event is Qo "- 
Multiplying by the probability that m neighbors are up, we finally get 



GAh)=y^ ] P*" go'"" HzJ - 2mJ - h). 
„ \m 



(3.11) 



m=0 



Chapter 3. Distribution of avalanche sizes on the Bethe lattice 



30 



Arguing similarly, we can write the equation for Gs{h) for s = 2, 3 etc. These equations 
simplify considerably when expressed in terms of the generating function G{x\h), and we 
get 

G(x\h) = xy f ^ 1 P*" QixY-'^ MzJ - 2mJ - h). (3.12) 

In numerical simulations, and experiments, it is much easier to measure the avalanche 
distribution integrated over the full hysteresis loop. To get the probability that an avalanche 
of size s will be initiated at any given site O in the interval when the external field is 
increased from hi to h2, we just have to integrate G{x\h) in this range. For any h, the 
value of dG/dx at x = 1 is proportional to the mean size of an avalanche, and thus to the 
average slope of the hysteresis loop at that h. 



3.2 Explicit calculation for the rectangular distribution 

While the general formalism described in the previous section can be used for any distri- 
bution, and any coordination number, to calculate the avalanche distributions explicitly, we 
have to choose some specific form for the probability distribution function. In this section, 
we shall consider the specific choice of a rectangular distribution: The quenched random 
field is uniformly distributed between — A and A, so that 

cf>{hi) = ^, for -A<h,<A. (3.13) 

In this case, the cumulative probabilities Pm{h) become piece wise linear functions of 
h, and /^-dependence of the distribution is easier to work out explicitly. We shall work out 
the distributions for the linear chain (z = 2), and the 3-coordinated Bethe lattice. 



3.2.1 The linear chain {z = 2) 

The simplest illustration is for a linear chain. In this case the self-consistent equation, for 
the probability P* [Eq. (2.9)] becomes a linear equation. This is easily solved, and explicit 
expressions for Qo, and Q{x) are obtained (see Appendix A.l). The different regimes 
showing different qualitative behavior of the hysteresis loops are shown in Fig. 3.2 

For h < 2 J — A (region A), all the spin remain down. For h > A, all spins are up 
(region D). For A < J, we get a rectangular loop and the magnetization jumps discontinu- 
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Figure 3.2: Behavior of RFIM in the magnetic field - disorder {h — A) plane for a linear chain. The 
regions A-D correspond to qualitatively different responses. In region A all spins are down and in 
region D all are up. The avalanches of finite size occur in region B and C. 



ously from — 1 to +lin a single infinite avalanche, and we directly go from region A to D 
as the field is increased. For A > J, we get nontrivial hysteresis loops. 

The hysteresis loops for different values of A = 0.5, 1.5 and 2.5 are shown in Fig. 3.3. 
If A is sufficiently large (A > J), we find that the mean magnetization is a precisely 
linear function of the external field for a range of values of the external field h (region B in 
Fig. 3.2). For larger h values, the magnetization shows saturation effects, and is no longer 
linear (region C). 

The explicit forms of the generating function Q{x) are given in the Appendix A.l. 
We find that in region B, the function Q{x) is independent of the applied field h. The 
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distribution function Gs{h) has a simple dependence on s of the form 

G,{h) = A^s (3.14) 
where Ai is a constant, that depends only on J/ A, and does not depend on s or h. 

In region C, the mean magnetization is a nonlinear function of h. But Q{x) is still a 
rational function of x. From the explicit functional form of Q{x) and G{x\h) are given in 
the appendix A.l, we find that Gs{h) is of the form 

Gs{h) = [A[s + A',](^^^ , for s>2. (3.16) 
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Here A[ and A'^ have no dependence on s but are explicit functions of h. 

Integrating over h from — oo to oo we get the integrated avalanche distribution D^, 

/oo 
G,{h)dh. (3.17) 
-oo 

It is easy to see from above that the integrated distribution Dg also has the form 

D, = [A2S + B2][£j , for 5>2, (3.18) 
where the explicit forms of the coefficients A2 and B2 are given in the Appendix A.l. 



3.2.2 The case 



The analysis for the case 2; = 3 is very similar to the linear case. In this case, the self- 
consistent equation for P*{h) [Eq. (2.9)] becomes a quadratic equation. The qualitative 
behavior of solution is very similar to the earlier case. Some details are given in Ap- 
pendix A.l. We again get regions A-D as before, but the boundaries are shifted a bit, and 
are shown in Fig. 3.4. As before, in region B, the average magnetization is a linear function 
of h, and the avalanche distribution is independent of h. 

We find that in regime B, the distribution of avalanche sizes is given by 



G,{h) = N 



(2s] 



-1)1(3 + 2) 

where is a normalization constant given by 



(3.19) 



N 



2A 



:i - J/ A) 



(J/A) 



(3.20) 



It is easy to see that for large s, Gs{h) varies as 



S 2 K 



(3.21) 



where 

K = 4{1- J/A)(J/A). (3.22) 

In region B, J/ A is always less than 1/3, and so this function always has an exponential 
decay for large s. 

In the region C, we find that the avalanche distribution is of the form 

(25)! 



Gs{h) = N' 



(3 _ 1)1(3 + 2)! 



(3.23) 
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Figure 3.4: Behavior of RFIM in the magnetic field - disorder Qi — A) plane for Bethe lattice of 
coordination number 3. The qualitative behavior in different regions A-D is similar to that of a 
linear chain (Fig. 3.2). 



where N' is a normalization constant independent of s, and k is a cubic polynomial in the 
external field h: 



1 



8(1 - 2J/A) 



[{9 - 53(J/A) + 119(J/A)2 - 107(J/A)3} 

+ {-5 + 10(J/A) + 11(J/A)2} (VA) 
+ {3-9(J/A)2}(/,/A)2 + (/,/A) 



(3.24) 



For any fixed s, the integrated distribution Dg can be evaluated explicitly, but become 
lengthy even for small s. 
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3.3 General distributions 

The analysis of the previous section can, in principle, be extended to higher coordination 
numbers, and other distributions of random fields. However, the self-consistent equations 
become cubic, or higher order polynomials. In principle, an explicit solution is possible 
for z < 5, but it is not very instructive. However, the qualitative behavior of solutions is 
easy to determine, and is the same for all z > i. We shall take z = Am the following for 
simplicity. Since we only study the general features of the self-consistent equations, we 









A \ 



hi 



Figure 3.5: A schematic plot of a unimodal random field distribution which asymptotically go to 

zero at iboo. 

need not pick a specific form for the continuous distributions of random field distribution 
(f>{hi). We shall only assume that it has a single maximum around zero and asymptotically 
go to zero at ±oo, as shown in Fig. 3.5. 

For small width (A) of the random field distribution i.e. for weak disorder the mag- 
netization shows a jump discontinuity as a function of the external uniform field, which 
disappears for a larger values of A (section 2.4). For fields h just lower than the value 
where the jump discontinuity occurs, the slope of the hysteresis curves is large, and tends 
to infinity as the field tends to the value at which the jump occurs. This indicates that large 
avalanches are more likely just before the first order jump in magnetization. 

For z = 4:, the self-consistent equation for P*{h) [Eq. (2.9)] is cubic 

aP*3 + 5p*2 ^ ^p* ^ J ^ (3.25) 
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where a, b, c and d are functions of the external field h, expressible in terms of the cumula- 
tive probabilities pi^i = to 3, 

a = P3 -^P2 + 3pi - Po, (3.26a) 
b = 3p2 - Qpi + 3po, (3.26b) 
c = 3pi — 3po — 1, (3.26c) 
d = po. (3.26d) 

This equation will have 1 or 3 real roots, which will vary with h. We have shown this 
variation for the real roots which lie between and 1 in Fig. 2.5 for the case where (f>{hi) 
is a simple distribution given by Eq. (2.1 1). 

We have also solved numerically the self-consistent equation for P* for other choices 
of p{hi), like the gaussian distribution, and for higher z{= 4, 5, 6). In each case we find that 
the qualitative behavior of the solution is very similar. Note that the rectangular distribution 
discussed in the previous section is very atypical in that both the coefficients a and b vanish 
for an entire range of values of h. 

In the generic case, we find two qualitatively different behaviors: For larger values of 
A, there is only one real root for any h. For A sufficiently small, we find a range of h 
where there are 3 real solutions. There is a critical value Ac of the width which separates 
these two behaviors. For the particular distribution chosen [Eq. (2.1 1)], Ac ~ 2.10382. 

In the first case, the real root is a continuous function of h, and correspondingly, the 
magnetization is a continuous function of h. This is the case corresponding to A = 2.5 in 
Fig. 2.5. 

For smaller A < Ac, for large ±h there is only one root , but in the intermediate 
region there are three roots. The typical variation is shown for A = 1.5 in Fig. 2.5. In the 
increasing field the probability P*{h) initially takes the smallest root. As h increases, at a 
value h = /^disc, the middle and the lower roots become equal and after that both disappear 
from the real plane. At h = hdisc the probability P*{h) jumps to the upper root. Thus 
for A < Ac there is a discontinuity in P*{h) which gives rise to a first order jump in the 
magnetization curve. 

The field /^disc where the discontinuity of magnetization occurs, is determined by the 
condition that for this value of h, the cubic equation [Eq. (3.25)] has two equal roots. The 
value of P* at this point, denoted by P^^c satisfies the equation 

3aoPf c + 2&oP;„c + Co = 0, (3.27) 
where ao, &o and cq are the values of a, 6 and cat h = /^disc- 
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We now determine the behavior of the avalanche generating function Gs{h) for large s 
and h near /^disc- The behavior for large s corresponds to x near 1. So we write x = 1 — S, 
with S small, and h = /^disc — £• Near hdisc, a^b^ . . . vary linearly with e and 

^*~n.c-«V^+0(e), (3.28) 

where a is a numerical constant. 

Since Q{x = 1) = 1 — P*{h) , if x differs slightly from unity Q{x) also differs from 
1 — P*{h) by a small amount. Substituting x = 1 — 5 and Q{x = 1 — S) = 1 — P* — F{e, 5) 
in the self-consistent equation for Q{x) [Eq. (3.9)], where both 5 and F axe small, using 
Eq. (3.27), we get to lowest order in 5, e and F 

+ ^^/^F - -f''5 = 0, (3.29) 

where (3 and 7 are some constants. Thus, to lowest orders in e and S, F is given by 

F = (l/2) 



^/32e + 472^ - [3y^e 



(3.30) 



Thus Q{x) has leading square root singularity at x = 1 + Consequently, G{x\h) 



47^ 

^2, 



will also show a square root singularity x = 1 + f^. This implies that the Taylor expansion 
coefficients Gs{h) vary as 

Gs{h) ~ s-l (1 + 0) . for large s. (3.31) 

At e = 0, we get 

a(/^disc) ~ 5-t. (3.32) 



Thus at h = /^disc the avalanche distribution has a power law tail. 

To calculate the integrated distribution D^, we have to integrate Eq. (3.31) over a range 

2 

of e values. For large s, only e < ^ contributes significantly to the integral, and thus we 
get 

_ 5 

Ds ^ s 2 ^ for large s. (3.33) 

Thus the integrated distribution shows a robust (—5/2) power law for a range of disorder 
strength A. 



Minor hysteresis loops 
on the Bethe lattice 



In this chapter, we derive exact self-consistent equations to obtain magnetization on the 
minor loops as a function of external field for arbitrary distribution of quenched random 
fields on a Bethe lattice. The return hysteresis loops for the linear chain was obtained by 
Shukla (2000). In sec. 2.4, we have discussed how to obtain the magnetization on the 
lower hysteresis curve, i.e. if we start with h = — oo, when all the spins are down, and 
then slowly increase the external field. Now suppose the system is on the lower hysteresis 
curve at some external field hi. Decreasing the field from hi to some field h2 and then 
again increasing to hi, we obtain the first minor loop. Similarly starting from the first 
minor loop at some field h^ and decreasing the field to h^, and then increasing to h^, we 
obtain the second minor loop and so on. Figure 4.1 shows two minor loops. In general, 
the n-th minor loop for n > 1 is obtained from the lower half of [n — l)-th minor loop by 
decreasing the field from /i2n-i to h2n and then increasing to /i2n+i < ^2n-i • This involves 
{hn} = hm hn-i, hi, the history of all the turning points from hi to hn- In the next 
section we will obtain the exact expressions for the magnetizations on the minor loops for 
arbitrary distributions of random fields. Similar results were later obtained independently 
by Shukla (2001). 



4.1 Magnetization on minor loops 

In sec. 2.4, we have determined the average magnetization in the deep inside the Cayley 
tree, on the lower hysteresis loop for arbitrary distributions of random field distributions. 
The average magnetization is equivalent to the magnetization at the root O of the Cayley 
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tree [Fig. 2.4], in the limit the number of generations, n — > oo. We obtained the magneti- 
zation at O as, 

Prob(3o = +l\h)=j2('] [P-'ihT [1 - P*!/^)]^"" Pmih), (2.10) 

m=0 \ / 

where P*{h) is the limiting value (r ^ n, and the limit n — > oo) of conditional probability 
P^'''\ that a spin on the r-th generation will be flipped when its parent spin at generation 
r — 1 is kept down, the external field is h, and each of its descendent spins has been relaxed, 
and is obtained by solving the polynomial equation 

p^h) = E " ^] [PHh)r [1 - PHh)r'-"^ p^{h), (2.9) 

m=0 \ / 

and Pm{h) is the probability that a spin flips up, given that exactly m neighbors are up, 
which is obtained by integrating the random field distribution (f>{hi) as, 

Pm{h) = / (pih,) dh,. (2.7) 

J(z-2m)J-h 



Similarly for the upper half of the hysteresis loop, when the external field is decreased 
from oo, we can define Q^'''\h) to be the conditional probability that a spin on the r-th 
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generation will be flipped down when its parent spin at generation r — 1 is kept up, the 
external field is decreased from oo to h, and each of its descendent spins has been relaxed. 
The limiting value Q*{h) also satisfies self-consistent equation 

Q'-ih) = E " ^] [1 - Q'-ihT [Q''{h)r'-'^ [I - P-r.+i{h)] , (4.1) 

m=0 \ / 

and in terms of Q*{h) the upper half of the major loop can be obtained. Since Pm+i {h) = 
Pm{h + 2J), the recursion relation satisfies by 1 — Q^'''\h — 2J), is same as the relation 
satisfies by P^'''\h) which is given by Eq. (2.8). Therefore, we conclude that g('')(/i-2 J) = 
1 - pW(/i). 

4.1.1 First minor hysteresis loop 

Suppose the system is on the lower hysteresis curve at some external field hi. Now if 
the field is decreased from hi to some field h2 and then again increased to hi, return point 
memory [section 2.1] ensures that the loop closes. This is the first minor loop [see Fig. 4.1]. 
Now when the applied field is increased from — oo to /^i and then decreased to a field /i2, to 
find out the spins which can flip down we need to consider only about the subset of spins 
which flipped up at field hi. Suppose a spin at a randomly chosen site flips up at field hi. 
As a result, the net local field at each of its nearest neighbors increases by an amount 2. J 
and some of down neighbors might become unstable. We flip up those spins at time step 
1 . After flipping them more of their neighbor might become unstable. We flip them up in 
time step 2 and so on. This process will be continued till the avalanche stop. Figure 4.2 
shows the order at which spins flip during a particular avalanche. Now in this avalanche 
if a spin Si flips up at time step t and as a result, if m of its neighbors flip at time step 
t + 1, then the local field at i will increase by 2m J. Therefore when the field is decreased 
to h2 > hi — 2.J, Si can not flip back at h2 unless all the neighbors which had flipped at 
time step t + 1 after Si flipped up, again flip back at h2. Therefore, the spin which was the 
initiator of the avalanche (which flipped at time step 0) can flip down at h2 only at the end, 
after all the spins of that avalanche flip back and in this flip-back avalanche the spins flip 
exactly in the reverse time order to the previous avalanche. This property will be called the 
time ordering property of spin-flip-back process. 

Consider the case, when the system is on the lower half of the major loop at field hi 
and then the field is decreased to hi — 2.J . Then all the neighbors of a vertex i which had 
flipped up at hi after Si flipped up will flip back at hi — 2.J, since they flipped up when 
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Figure 4.2: Time order at which spins flip during a particular avalanche. 



their local fields had been increased by 2 J. Therefore, the conditional probability that a 
spin is down at hi — 2J, given its parent spin is up is same as the conditional probability 
that a spin is down at hi, given its parent spin is down, which is 1 — P^'''\hi). The later 
is again equal to Q^'''\hi — 2J), the conditional probability that a spin is down when the 
field is decreased from oo to /^i — 2 J, given its parent spin is kept up. Therefore the reverse 
magnetization curve starting from hi, meets the upper major half at hi — 2 J and merge with 
it for h2 < hi — 2J. This result can be generalized for arbitrary graphs, which is discussed 
in section 4.2. Thus, we can consider the first minor loop in the range [hi — 2 J, hi]. Since 
in this range of external field the spin-flip-back process obeys time ordering, if a spin Si 
flips up at hi and flips down at h2, then the probability that a neighbor of it at generation 
r is up before Si flips back at h2 is same as the neighbor was up before Si flipped up at hi, 
given by P^'''\hi). The probability that a neighbor is down before Si flips down at Si can 
be splitted into two parts: 

(1) it didn't flip up after Si flipped up at /^i and 

(2) it flipped up after Si flipped up at hi and flips back at h2 before Si flips back. 
Consider a site X at some level r of the Cay ley tree [Fig. 3.1]. We call the subtree 

formed by X and its descendents Tx, the subtree rooted at X. We keep its parent spin Y 
at generation r — 1 down, and relax all the sites in Tx at the uniform field hi. Let R^^\hi) 
be the probability that sx remains down after sy turned up at hi. For r ^ n, in the limit 
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n — > oo, these probabilities tends to limiting value R*L{hi), given by 

R'.ih,) = 1 - E " [P'-ih.T [1 - P^l/^l)]^-'-" P-ra + l{h,). (4.2) 

— n \ m 

m=0 \ / 

Let G^^\h2, hi) be the conditional probability that: 

(a) sx was down at hi, given sy was down, 

(b) Sx flipped up at hi after sy flipped up and 

(c) Sx flips back at h2, given sy is still up. 

Then a recursion relation for G^^\h2, hi) in terms of G^^~^^\h2, hi) can be written as 

G^:\h2.hi) = j: ('-^) [p^{hi)r [R^:^'\hi) + G^:^'\h2.hi)Y-'-"' 

X [p^+i{hi) - p^+i{h2)] , (4.3) 
and its limiting value G* (/i2, ^i), satisfies the self-consistent equation 

GUh2,hi) = 2 ('~^) [P*(/^i)r [RUhi) + GUh2,hi)Y~'~"' 

X Ipm + lihi) ~ Pm + 1 ih2)] . (4.4) 

This is a polynomial equation in G* (/i2, /^i) of degree z — 1, whose coefficients are func- 
tions of hi and h2 through P*{hi), R*_{hi), Pm{hi) and p.m{h2). To determine G* (/i2, ^i) 
for any given pair of external fields hi and h2, we have to first solve the self-consistent 
equation for P*{hi) [Eq. (2.9)]. This then determines R^{hi) using Eq. (4.2), and then, 
given P*{hi) and R^{hi), we solve for G* (/i2, hi) by solving the [z — l)-th degree poly- 
nomial equation Eq. (4.4). Now the decrease in magnetization, when the field is decreased 
from hi to h2, is determined by the probability that a spin at O was up at hi and turns down 
at h2, given by, 

Prob(so = -1; h2\ So = +1; hi) = 

E (')[P''ihi)r[RUhi) + GUh2,hi)Y~'^ [Pn^ihi)-p^ih2)]. (4.5) 

m=0 y^J 

This determines the upper half of first minor loop. 

Similarly when the field is again reversed from h2 to h^ < hi, using the symmetry 
between up an down spins it is easy to see that again the time ordering property holds. 
Therefore the probability that the neighbor of a spin Si is down before Si flips up at h^ is 
[R^^\hi) + G*^^(/i2, hi)]. The probability that a neighbor is up at h^ before Si flips up is 
given by sum of two probabilities i?+^(/i2, hi) and G^^\hs, h2,hi); where i?+^(/i2, hi) is 
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the probability that: (a) sx is up at hi, given that sy is kept down and Tx is relaxed, (b) 
sy flips up at hi and (c) sx remains up after sy flips down at h2 and Tx is relaxed and 
G^l\h^^ h2,hi) is the probability that: (a) sx is up at hi, given that sy is kept down and 
Tx is relaxed, (b) sy flipped up at hi, (c) sx flipped down after sy flipped down at h2 and 
(d) Sx flips back at h^, given sy is still down. 

i?+ ^ (/i2 , /ii ) is the equal to the probability that the spin is up at /^i minus the probability 
that it becomes down at h2. Its limiting value is given by 

RX{h2M) = P'-ih) - E [P'^ihiT [R-ihi) + G'L{h2,hi)] 

m=0 \ ITT- J 

X IPmihi) - Pmih2)] . (4.6) 
The limiting value G^{h3^ /i2, ^i) satisfies the self-consistent equation 

Gl{hs,h2,hi)= J2 [Rl{h2,hi) + Gl{hs,h2,hi)]'^ 

m=0 \ ITT- J 

R^hi) + (/.2,/^l)] [Pm{h:,) - p^{h2)] . (4.7) 

Solving the above self-consistent equation [Eq. (4.7)] we determine (^^(/is, /i2, /^i) and 
then the increase in magnetization, when the field is increased from h2 to h^, is determined 
in terms of the following probability: 

Prob(so = +1; h^\so = -1; h2) = 

z / \ 

{ Z \ V ~\ TTh r ~\ Z — TTh 

E ][RX{h2,hi) + GX{h^,h2,hi)\ [i?*(/.i) + G*(/.2,/^l)] 

m=o y^J 

xlPmih^) - Pm{h2)], (4.8) 

which determines the lower half of first minor loop. 



4.1.2 General minor hysteresis loops 

In the previous sub section, we obtained the first minor loop. The other minor loops can be 
obtained similarly. In all the minor loops the spin-flip-back process obeys time ordering. 
In general, the n-th minor loop for n > 1 is obtained from the lower half of (n — l)-th 
minor loop by decreasing the field from /i2n-i to h2n and then increasing to /i2n+i < ^2n-i • 
For convenience, we will use the notation {hn} = h.^ hn-i , /^i for the history of all the 
turning points from hi to hn- 

On the upper half of the n-th minor loop (n > 1), when the field is decreased from 
h2n-i to h2n, thc probability that a neighbor of a spin is up before (which is deep 
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inside the tree) flips down at is [R^{{h2n-2}) + G'^({/i2n-i})]- The probability that a 
neighbor of Si is down before it flips down is given by [i?* ({/i2n-i}) + G* ({/i2n})], where 
E*_{{h2n-i}) is given by, 



- E " ^] [Ki{h2n-2}) + Gl{{h2n-,})r 

m=0 \ ITT- / 

x[RU,{h2n-3}) + GU,{h2n-2})r-'-"' 

x[Pm + l{h2n-l) — Pm + 1 ih2n-2)], (4.9) 

and G* ({/i2n}) satisfies the self consistent equation 

GU{h2n}) = J: (" ~ %Rl{{h2n-2}) + Gl{{h2n-,})r 

x[R^_{{h2n-l}) + G^_{{h2n})r-'-''^ 

ih2n)]. (4.10) 

Therefore the decrease in magnetization, when the field is decreased from /i2n-i to h2n, is 
obtained from 

Prob(so = -l;h2n\ so = +l;/i2n-i) 

= E ("][Rl{{h2n-2}) + Gl{{h2n-l})r 
m=0 V^J 

x[RU{h2n-l}) + G*_{{h2n})Y-''^ K (/^2.-l ) " (/^2.)] . (4.11) 

Similarly on the lower half, the increase in magnetization, when the field is increased 
from h2n to /i2n+i, is obtained from 

Prob(so = +I;h2n+i \ so = -l;h2n) 
I m 



= E m{{h2n}) + Gl{{h2n+l})y 

m=0 V^J 

x[RU{h2n-l}) + G*_{{h2n})Y-''^ K (/^2.+l ) " (/^2.)] , (4.12) 

where R\{{h2n}) is given by 

Rl{{h2n}) = [Rl{{h2n-2}) + Gl{{h2n-l})] 

- E " ^] [Ki{h2n-2}) + Gl{{h2n-l})r 

m=0 \ ITT- / 

x[R*_i{h2n-l}) + G*_i{h2n})r-'-"' 

X[pm{h2n-l) - Pm{h2n)], (4.13) 
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and G^{{h2n+i}) is obtained by solving the self consistent equation 



Gl{{h2n+l}) 



z-1 



E 



( 



z - 1 



m 



) 



[Rl{{h2n}) + Gl{{h,^+^})] 



x[RU{h2n-l}) + GU,{h2n})] 
x[Pm{h2n+l) - Pm{h2n)]- 



Z — l — m 



(4.14) 



In Fig. 4.1, we have plotted first two minor hysteresis loops, generated by solving 
above equations, for three coordinated Bethe lattice. The random field distribution is given 
by Eq. (2.11) and we choose A = 1.5. 



4.2 Merging of different stable configurations 

In the previous section, we see that the two ends of minor loop are at major loop, with the 
external field is differed by 2J. In this section we generalize this to any stable configuration 
on any graphs. We prove that, for RFIM on a connected general graph G, for a given 
realization of random fields, all the stable configurations at external field h go to unique 
stable configuration C{h ± 2z*J) when the field is monotonically changed to h ± 2z*J; 
where z* is the minimum number required such that, any connected subgraph g of G has 
at least one vertex such that, the number of edges in g connected to that vertex is < z*. For 
example, z* = 2 for square lattice, and for Bethe lattice z* = 1. 

Proof. — Consider two stable configurations Ci{h) and C2{h) at field h. We can de- 
compose vertices of the graph G, into sets: (1) A„„, up-spins in both configurations, (2) 
Aud, up-spins in Ci and down-spins in C2, (3) Adu, down-spins in Ci and up-spins in C2 
and (4) Add, down-spins in both configurations. Consider when the external field is in- 
creased monotonically by 2z*J. Since in the zero temperature dynamics, in the increasing 
field field the spins flip only once, and the order in which various spins are relaxed does 
not matter, we can increase the field in one step to h + 2z*J and then first relax spins from 

sets Aud and Adu- Now the set Aud can be written as union of disjoint subsets A^^], A^^] 

Consider one such subset, which is on a subgraph g. On this subgraph, the local field at a 
vertex i in configurations Ci and C2 at field h are. 



liCrih)) = ^9 J _ ^'^J + > 0, 



(4.15) 



(4.16) 
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where zf are the number of vertices in g connected to vertex i, z\ are the number of vertices 
in the set A^^ connected to i and ji is the contribution to the local field from the sets A„„ 
and Au. Since £f '^''^^ - 4^'^''^^ > 0, from Eq. (4.15) and Eq. (4.16) we get - ^' > 
or zl > 1 as 2;' > 0. This means that, there can not be a subset of A„d which has only one 
element. 

Now when the external field is increased to /i + 2z^ J , the local field at vertex i in 
configuration C2 becomes, 

<• = {K + + z^J - z\J + /,} + {(4 - zl)J + 2z\J} . (4.17) 

Therefore, at all those vertices in g where zf < z* (there is at least one such vertex in g, 
by the definition of z*), the local field will become positive. So, the spins in C2, at those 
vertices will flip up and the original subset will will shrink to a new one on a different 
reduced subgraph g, and the same argument holds for it also. Therefore, after iterative use 
of this relaxation procedure, all the subsets of Aud will become null sets and by the same 
argument, it is also true for Adu- Therefore, after relaxing all the spins from the sets Aud 
and Adu, the resultant unstable configurations Ci{h + 2z*J) and C2{h + 2z*J) are identical, 
and hence relaxation of the remaining unstable spins from the set Add, will lead to the same 
final stable configuration. From the symmetry between up and down spins, it is obvious 
that the configuration Ci{h) and C2{h) go to the same final configuration, when the field 
is decreased by 2 2;* J. 

For a square lattice z* = 2, as in any connected subgraph of it, there exists at least 
one vertex form which the number of edges connected to the subgraph is < 2. Therefore 
any two different stable configurations should merge to one configuration, when the field 
is increased or decreased by 4J. In Fig. 4.3, we consider two stable spin configurations 
(a) and (b) at external field h = 0, on a square lattice, for a given realization of random 
fields. The lattice size is 50 x 50. The spins which are down in (a) and up in (b) are shown 
in (c) in black color. Now when we increase the external field to 4J, configurations (a) 
and (b) evolve to new stable configurations (d) and (e) respectively. We see that these two 
configurations (d) and (e) are identical as seen from their difference configuration (c). 



Chapter 4. Minor hysteresis loops on the Bethe lattice 



47 




(d) (e) (f) 

Figure 4.3: (a) and (b)are two different stable spin configurations at /i = 0, witli tlie same realiza- 
tion of random fields, on a square lattice of size 50 X 50. The up spins ai^e represented by black and 
down spins by white color, (c) shows the difference between (a) and (b). The spins, which are down 
in (a) and up in (b) are represented by black color in (c). (d) and (e) are new stable configurations 
obtained from (a) and (b) respectively, ath = 4 J. (f) shows the difference between (d) and (e). 



Chapter 5 



Hysteresis on regular lattices 
in the low disorder limit 



In general, Bethe approximation is expected to work well for noncritical properties. Is the 
Bethe approximation is a good approximation for regular lattices? This is the question we 
address in this chapter. Surprisingly, for asymmetrical distribution, the answer can be no. 

In this chapter, we discuss the low disorder limit of the hysteresis loop in the random 
field Ising model (RFIM) on periodic lattices in two and three dimensions. We find that 
the behavior of hysteresis loops depends nontrivially on the coordination number z (Sab- 
hapandit et al. 2002). For 2; = 3, for continuous unbounded distributions of random fields, 
the hysteresis loops show no jump discontinuity of magnetization even in the limit of small 
disorder, but for higher z they do. This is exactly as found in the exact solution on the 
Bethe lattice (Dhar et al. 1997). 

As discussed in the introduction, random field Ising model was first studied in the 
context of possible destruction of long range order by arbitrarily weak quenched disorder 
in equilibrium systems. Accordingly the distribution of random field was assumed to be 
symmetrical. However, in hysteresis problem, the symmetry between up and down spins 
state is already broken by the specially prepared initial state (all down in our case), and the 
symmetry of the distribution plays no special role. 

The analytical treatment of self-consistent equations on the Bethe lattice is immediately 
generalized to asymmetrical case. However, we find that for asymmetrical distributions 
the behavior of hysteresis loops in euclidean lattices can be quite different from that on 
the Bethe lattice. On hypercubical lattices in d dimensions, there is an instability related 
to bootstrap percolation, that is absent on the Bethe lattice. This reduces the value of 
the coercive field /icoer away from the Bethe lattice value 0{ J) to zero, where J is the 
exchange coupling. We note that the limit A — > is somewhat subtle, as the system size 
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L* required for self-averaging diverges very fast for small A, and the finite-size corrections 
to the thermodynamic limit tend to zero very slowly. 

In the following, we shall assume that the distribution has a asymmetrical shape, given 

by 

^{hi) = ^expi-h,/A)9ihiy, (5.1) 

where 9 is the step function. The mean value of hi can be made zero by a shift in the 
value of the external uniform field. Our treatment is easily extended to other continuous 
unimodal distributions. The exact form of (f>{x) is not important, and other forms like 
exp(— x — e~^) which fall sharply for negative x have the same behavior. 

For a given distribution (f>{hi), we define Pm{h) with < m < 2; as the conditional 
probability that the local field at any site i will be large enough so that it will flip up, if m 
of its neighbors are up, when the uniform external field is h. Clearly 

Pm{h) = / (t){hi) dh,. (2.7) 

J (z-2m)J-h 

Clearly, for any given value of h, the magnetization depends on the distribution (j){hi) 
only through Pm{h). 



5.1 Hysteresis on three coordinated lattices 

Consider first the case of the two-dimensional hexagonal lattice with z = 3. For periodic 
boundary conditions, if A = 0, starting with a configuration with all spins down, clearly 
one has hcoer = 3 J. For A 7^ 0, the site with the largest local field flips first, and then if 
h > J, pi{h) = 1, this causes neighbors of the flipped spin to flip, and their neighbors, 
and so on. Thus, so long as there is at least one flipped spin, all other spins also flip, and 
the magnetization is 1. The largest local field in a system of L'^ spins is of order 2 A In L. 
Once this spin turns up, other spin will flip also up, causing a jump in magnetization from 
a value ~ —1 to a value +1 in each sample. Hence the coercive field, (the value of h where 
magnetization changes sign) to lowest order in A, is given by 

h^^^, = 3 J - 2A In L, for 1 < In K J/ A. (5.2) 

Sample to sample fluctuations in the position of the jump are of order A. On averaging 
over disorder, the magnetization will become a smooth function of h, with the width of the 
transition region being of order A. 
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{h-J)/A 

Figure 5.1: Magnetization in the increasing field. The curves for the two values of A coincide. 
Curves A is for hexagonal lattice of size 4096^ and B is for a three coordinated lattice in three 
dimensions [see Fig. 5.2] of size 256"^. 

For a fixed A ^ J, if L is increased to a value near exp(J/A) = L^^^, h^oer decreases 
to a value near J. For h ^ J, pi{h) is no longer nearly 1, but po{h) ~ 0, and p2{h) ~ 
P3{h) ~ 1. The value of magnetization depends only on pi{h), which is a function of 
h = [h — J) / A. As h increased from — oo, pi{h) increases continuously from to 1. 

Note that for A = 0.0 IJ, L^^-^ ~ lO"*"^. Therefore it is impossible to study the large L 
limit with the available computers. To avoid the problem of probability of nucleation being 
very small for h near J, we made the local field at a small fraction of randomly chosen sites 
very large, so that these spins are up at any h. The number of such spins we choose to be 
of order L, so that their effect on the average magnetization is negligible. Introduction of 
these "nucleation centers" makes L* 0{\/Tj) ( the average separation between centers), 
and hcoer drops to a value near J, so that, we can study the large L limit with available 
computers. For L > L^^^, the behavior of hysteresis loops becomes independent of L. 

In Fig. 5.1, curve A shows the result of a simulation on the hexagonal lattice with 
L = 4096, and periodic boundary condition. We see that magnetization no longer under- 
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goes a single large jump, but many small jumps. In the figure, we also show the plot of 
magnetization when the random field at each site is decreased by a factor 10. This changes 
the value A from O.lJtoO.OlJ. However, plotted as a function of h, the magnetization for 
these two different values (for small A) fall on top of each other /or the same realization 
of disorder (except for the overall scale A). Thus we can decrease A further to arbitrarily 
small values, and the limit of A is straightforward for each realization of disorder. 
Then, averaging over disorder, for a fixed A, we see that /^coer tends to the value J as A 
tends to 0. Also, we see that there is no macroscopic jump-discontinuity for any non-zero 
A. 





Figure 5.2: A three coordinated lattice {z = 3) in three dimensions. 



We also show in Fig. 5.1 [curve B], the results of simulation of a 3-dimensional lattice 
with z = ?> [shown in Fig. 5.2] of size 256"^ with periodic boundary condition. The behavior 
is qualitatively same as that in two dimensions. The value of /^coer = in the limit A — > 
is same for symmetrical distribution, and also is the same as predicted by the Bethe 
approximation. 
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5.2 Bootstrap instability in RFIM on square lattice 

On the square lattice also, the value of /^coer is determined by the need to create a nucleation 
event. Arguing as before, we see that /^coer to lowest order in A is given by 



hcoer ~ 4J — 2Alniv, 



for l<lniv<J/A. 



(5.3) 



Adding a small number of nucleation sites suppresses this slow transient, and lowers hcoer 
from 4 J to a value near 2J. However, in this case, even after adding the nucleation centers, 
the system shows a large single jump in magnetization, indicating the existence of another 
instability. We observed in the simulation that at low A, as h is increased, the domains of 
up spins grow in rectangular clusters [see Fig. 5.3] and at a critical value of /^coer, one of 
them suddenly fills the entire lattice. This value hcoer fluctuates a bit from sample to sample. 



Figure 5.3: A snapshot of the up-spins just before the jump (h = 1.998243 J). The lattice size is 
200 X 200 and A = O.OOIJ. Initial configuration is prepared with 0.05% up-spins. 
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Figure 5.4: Distribution of the scaled coercive field on a square lattice for different lattice size L^. 



In Fig. 5.4 we have plotted the distribution of the scaled variable he = (/^coer — 2 J)/ A for 
different system sizes L, for A = O.OOIJ. The number of different realizations varies from 
10'* (for the largest L) to 10^ (for the smallest L). Note that the distribution shifts to the left 
with the increasing system size, and becomes narrower. 

This instability can be understood in terms of bootstrap percolation process BPm (see 
Adler 1991, for a review). Bootstrap percolation was first considered by Chalupa et al. 
(1981) (also Kogut and Leath 1981) and was subsequently studied by many others in a 
variety of contexts. The process BPm is define as follows: On a (/—dimensional lattice, 
sites are independently occupied with a probability p and the resulting configuration is 
taken as the initial configuration, which is evolved by the following rules: 

(a) the occupied sites remain occupied forever, 

(b) an unoccupied site having at least m occupied neighbors, becomes occupied. 

For m = 2, on a square lattice, in the final configuration, the sites which are occupied 
form disjoint rectangles, like the cluster of up-spins in Fig. 5.3. It has been proved that 
in the thermodynamic limit of large L, for any initial concentration p > 0, in the final 
configuration all sites are occupied with probability 1 (Aizenman and Lebowitz 1988). 

In the random field Ising model on a square lattice, for the asymmetric distribution 
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[Eq. 5.1] for h > 0, pm = I for m > 2, and any spins with more than one up-neighbors 
flips up. Therefore, stable clusters of up spins are rectangular in shape. The growth of 
domains of up spins is same as in the bootstrap percolation process BP2. 

Consider a rectangular cluster of up spins, of length / and width m. Let P(/, m) be 
the probability that, if this rectangle is put in a randomly prepared background of density 
Pi{h), this rectangle will grow by the BP2 process to fill the entire space. The probability 
that the random fields at any sites neighboring this rectangle will be large enough to cause 
it to flip up is pi{h). The probability that there is at least one such site along each of two 
adjacent sides of length / and m of the rectangle is (1 — g')(l — g"), where q = 1 — pi{h). 
Once these spins flip up, this induces all the other spins along the boundary side to flip up 
and the size of the rectangle grows to(/ + l) x (m + 1). Therefore 

P(/, m) > (1 - g')(l - g'")P(/ + 1, m + 1). (5.4) 

Thus the probability of occurrence of a nucleation which finally grows to fill the entire 
lattice is 

00 

Pnuc > P^{h)X{{l-q'f 

3 = 1 

00 

~ Po(/i) n [1 - exp{-pi(/i)j}]^ 

po(/i)exp f --^^ ) , forsmaU pi(/i). (5.5) 

The condition to determine /^coer is that for this value of h, P^^c becomes of order l/P^, 
so that we get 

Po(/ieoer)exp ( -— r)~7^- (5.6) 

This equation can be solved for /^coer for any given L. For the distribution given by Eq. (5.1), 
this becomes 

StT / ^coer H~ 2</ 



/ ^coer 4</ \ 

exp I 1 exp 



3-'"n — A 



for /ieoer < 2 J. (5.7) 



Therefore, the leading L-dependence of /^coer, to lowest order in A is given by 

3 



hnr\pr ~ '^J 111 



,(lnL- J/A) 



for J/A < In K exp(2J/A). (5.8) 



This agrees with our observation that the scaled critical field he shifts to the left with 
increasing system size. The width of the distribution of over which the coercive field 
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varies, can be calculated from the width over which the probability of having at least one 
nucleation in the entire lattice, i.e. 1 — (1 — Pnuc)^ ~ 1 — exp ( — PnucL^), changes from 
almost zero to almost unity: 

SK,,r ~ 7^ (5.9) 

in L 

Therefore, for any fixed A > 0, the jump will smeared out on averaging over disorder. Only 
in the limit A — > and L — > oo, the average magnetization will show a jump discontinuity. 

To test the validity of Eq. (5.6) in simulations, we put po{h) = 0.005 independent of h. 
Eq. (5.6) then simplifies to 

Pl(/icoer) ~ 777— T- (5-10) 

D in L 

In Fig. 5.5, we have plotted pi for the mean hcoer from Fig. 5.4 versus 1/ In L. The graph 
is approximately a straight line, which agrees with Eq. (5.10). The slope of the line is 
0.765 ± 0.009, less than in Eq. (5.10), which only gives an upper bound to /^coer- 

If h > 0, we will have P2 = I, and bootstrapping ensures that so long as po > 0, we 
will have all spins up in the limit of large L. This implies that hcoer = in this limit. 

If there are sites with large negative quenched fields, the bootstrap growth stops at 
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such sites. Hence the bootstrap instability cannot be seen for symmetric distributions. 
For a symmetrical distribution of random fields, the average distance from a nucleus to the 



nearest spin, which does not flips up even if it has two up neighbors, is L* = l/yl — p2{h). 
Therefore the average area covered by a nucleus is 

L 

L''l[{l-q'Y for L<L* (5.11a) 

and 

L*^l[{l-q'Y for L>L\ (5.11b) 

The condition to determine hcoer is that for this value of h, the average area due to the 
growth becomes 0{L'^), i.e. 

Po{hcoer)L^ X average area covered by a nucleus ~ 0{L'^). (5.12) 

From Eq. (5.1 lb) and Eq. (5.12), in the limit L — > oo, 

/fr\ na-g^T^^a)- (5-13) 

Nowpo(2J) = 1 -p2(2J),and the product nj=i(l - q^ f is 0(1) at h = 2J aspi(2J) = 
1/2. Therefore, /^coer ~ 0{2J). 

Even if the quenched fields are only positive, the instability does not occur on lattices 
with z = 3. On such lattices, if the unoccupied sites percolate, there are infinitely extended 
lines of unoccupied sites in the lattice. These cannot not become occupied by bootstrapping 
under BP2. Thus the critical threshold for BP2 on such lattices is not 0. 



5.3 Bootstrap instability in RFIM on cubic lattice 

The arguments for large void instability can be easily extended to higher dimensions. In 
d = 3, if h > 0, then Pm{h) = 1 for m > 3, therefore the spin flip process is similar to the 
spanning process of three dimensional BP3 (Cerf and Cirillo 1999). In this case, it is known 
that for any initial non-zero density, in the thermodynamical limit, the final configuration 
has all sites occupied with probability 1. The clusters of up-spins grow as cuboids, and at 
each surface of the cluster, the nucleation process is similar to that in two dimension. Let 
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e be the probability that, a nucleation occurs at a given point of a surface of the clusters of 
up spins which sweeps the entire two dimensional plane at h. 

e^Mh)e.p^-^y (5.14) 

The probability that, there exist at least one nucleation which sweeps the entire plane of 
size / X /, is 1 — (1 — e)'^ Therefore, the probability Pnuc, that a nucleation sweeps the 
entire three dimensional lattice at h satisfies 



Pnuc > Po{h)l[[l-{l-ef 

1=1 

oo 

~ Po{h) n [l - exp(-e/ 



3 



/=1 

exp(— A/a/e), for small e; (5.15) 
where A = ^^/^T({3/2). hcoer is determined by the condition that Pnuc niust be of the order 

A ( TT 



Po(/icoer)exp 



. ; exp . 

V6P2(/^coer) 

The leading P-dependence of hcoer is different in different ranges of /icoer, depend- 
ing on whether the strongest dependence of the left-hand side comes from variation of 
Po{h),pi{h) orp2{h). 

In the range 4 J < hcoer < 6 J: pm = I, for m > 1. Then we must have po(^coer) ~ 
1/P"^, which for the distribution given by Eq. (5.1) results 

/icoer ~ 6J - 3AlnP. (5.17a) 

The corresponding range of P, for the validity of of above equation is 1 ^ In P ^ 

(2J/3A). 

In the range 2J < hcoer < 4J: Pm = 1, for m > 2. Then in Eq. (5.16) the left hand side 
varies as exp [—A'/ ^Jpi{ hcoer)] , which gives 

/icoer~4J-2AlnMnP-— j , (5.17b) 

which is valid in the range (2 J/3 A) < In P < exp(2 J/A). 

In the range < hcoer < 2J: = 1 for m > 3. Then from Eq. (5.16), to the lowest 
order in A, we get 

/icoer ~ 2J - Alnln (^InP - , (5.17c) 



^ Ts- (5-16) 
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forexp(2J/A) < InK exp(exp(2J/A)). 

In the limit L ^ L*^^ = exp(exp(exp(2J/A))), the loop becomes independent of L, 
with hcoer 0. We have also verified the existence of jump in numerical simulation for 
z = 4: (diamond lattice) in three dimensions. 



Chapter 6 



Discussion 



Analytical treatment of problems having quenched disorder is usually difficult. There are 
few models having nontrivial quenched disorder that can be solved exactly. In this thesis, 
we set up exact self-consistent equations for the avalanche distribution function for the 
RFIM on a Bethe lattice. We were able to solve these equations explicitly for the rectan- 
gular distribution of the quenched field, for the linear chain z = 2, and the 3-coordinated 
Bethe lattice. For more general coordination numbers, and general continuous distributions 
of random fields, we argued that for very large disorder, the avalanche distribution is ex- 
ponentially damped, but for small disorder, generically, one gets a jump in magnetization, 
accompanied by a square-root singularity. For field-strengths just below corresponding to 
the jump discontinuity, we showed that the avalanche distribution function has a power-law 
tail of the form s~^^'^. The integrated avalanche distribution then varies as for large 
s. 

We have also studied the behavior the return loop, when the external field is increased 
from — oo to some value hi, and then decreased to a lower value h2 and again increased 
to the previous extremum value hi. We set up exact self-consistent equations to determine 
the magnetizations on all minor loops for arbitrary distributions of random fields. 

Some unexpected features of the solution deserve mention. Firstly, we find that the 
behavior of the self-consistent equations for 2; = 3 is qualitatively different from that for 
z > 3. The behavior for the linear chain (z = 2) is, of course, expected to be different 
from higher z. One usually finds same behavior for all 2; > 2. Mathematically, the reason 
for this unusual dependence is that the mechanism of two real solutions of the polynomial 
equation merging, and both becoming unphysical (complex) is not available for z = 3. 
Here the self-consistency equation is a quadratic, and from physical arguments, at least 
one of the roots must be real. That a Bethe lattice may show non-generic behavior for 
low coordination numbers has been noted earlier by Ananikyan et al. (1994) in their study 
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of the Blume-Emery-Griffiths model on a Bethe lattice. These authors observed that the 
qualitative behavior for 2; < 6 is different from that for 2; > 6. To find out whether this 
unusual z dependence of the hysteresis loop persists on regular lattices, we study it on the 
regular lattices in two and three dimension, in the limit of low disorder. We find that for 
asymmetrical distributions of random fields, there is a instability which is not present in 
three coordinated lattices, and hence the hysteresis curve is continuous for such lattices. 

The second point we want to emphasize is that here we find that the power-law tail in 
the distribution function is accompanied by the first-order jump in magnetization. Usu- 
ally, one thinks of critical behavior and first-order transitions as mutually exclusive, as 
first-order jump pre-empts a build-up of long-ranged correlations, and all correlations re- 
main finite-ranged across a first-order transition. This is clearly not the case here. In fact, 
the power-law tail in the avalanche distribution disappears, when the jump disappears. A 
similar situation occurs in equilibrium statistical mechanics in the case of a Heisenberg 
ferromagnet below the critical temperature. As the external field h is varied across zero, 
the magnetization shows a jump discontinuity, but in addition has a cusp singularity for 
small fields^ But in this case the power-law tail is seen on both sides of the transition. 

Note that for most values of disorder, and the external field, the avalanche distribution 
is exponentially damped. We get robust power law tails in the distribution, only if we 
integrate the distribution over the hysteresis cycle across the magnetization jump. But, in 
this case, the control parameter h is swept across a range of values, in particular across 
a (non-equilibrium) phase transition point! In this sense, while no explicit fine-tuning 
is involved in an experimental setup, this is not a self-organized critical system in the 
usual sense of the word. Recently Pazmandi et al. (1999) have argued that the hysteretic 
response of the Sherrington-Kirkpatrick model to external fields at zero temperature shows 
self-organized criticality for all values of the field. However, this seems to be because of 
the presence of infinite-ranged interactions in that model. 

In chapter 3, we discussed the behavior of avalanche distribution for various distribu- 
tions of random fields. A general question concerns the behavior of the avalanches for 
more general probability distributions. Clearly, if p{hi) has a discrete part, it would give 
rise to jumps in pi as a function of h, and hence give rise to several jumps in the hysteresis 
loop. These could preempt the cusp singularity mechanism which is responsible for the 

t Below Tc, the magnetization goes as, 

mr~.sign{h)\mo + A\h\^''-^'>^^], as h^O, 2<rf<4 (see Parisi 1988). 
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power-law tails. If the distribution is continuous, but multimodal, then it is possible 
to have more than one first order jump in the magnetization^ This is confirmed by ex- 
plicit calculation in some simple cases. If p{hi) has power-law singularities, these would 
also lead to power-law singularities in pi, and hence in P*{h). Even for purely continuous 
distributions, the merging of two roots as the magnetic field varies need not always occur. 
For example, it is easy to check that for the rectangular distribution, even for 2; > 4, we 
do not get a power law tail for any value of A. The precise conditions necessary for the 
occurrence of the power-law tail needs to be investigated further. 

Finally, we would like to mention some other open questions. Our analysis relied heav- 
ily on the fact that initial state was all spins down. Of course, we can start with other initial 
conditions. For example, start with the equilibrium state at temperature To and field ho, and 
then quench to zero temperature. Our present treatment cannot be applied to these cases 
as finite temperature brings about very nontrivial coorrelations between spins. It would be 
interesting to set up self-consistent field equations for them. In case of minor loops also, 
we have always started with h = —00 and then vary the field cyclically. Moreover, to 
find the magnetization at some particular point of the hysteresis curve, we start with the 
previous extremal field and change the field to the new value in one jump and then relax 
the system. It would be useful to find out some dynamical relations by which system can 
be evolved from any state by changing the field in infinitesimal steps. 

Another extension would be to make the rate of field-sweep comparable to the single- 
spin flip rate (still assuming zero temperature dynamics). This would mean some large 
avalanches in different parts of the sample could be evolving simultaneously. Then one 
could study the sweep-rate dependence of the hysteresis loops, and the frequency depen- 
dence of the Barkhausen noise spectra. This is perhaps of some relevance in real experi- 
mental data, and would also make contact with other treatments of Barkhausen noise that 
focus on the domain wall motion. 

Another case of some interest is other type of disorder e.g. the site-dilution case dis- 
cussed by Tadic (1996). It seems plausible from the structural stability of the mechanism 
which leads to the cusp singularity just before the jump-discontinuity in magnetization, 
in our model, introduction of site dilution would not change the qualitative behavior of 
solutions. 

We hope that many of these issues will be resolved in the next few years. 

^This would happen if P*{h) as a function of h shows a 'double S' curve. Then there must be at least 4 
values of h for which the slope of the curve is infinite. This is possible only if the equation determining P^^^^ 
[variant of Eq. (3.27)] is at least a quartic, hence only if z > 6. 
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A.1 Avalanche distribution on a linear chain 

For the case of a linear chain, the self-consistent equation, for the probability P* [Eq. (2.9)] 
is a linear equation, whose solution is, 

P^h) = T--^ -y (A.l) 

[1 - {pi - po)\ 

For h < 2J — A, po is zero, and hence P*{h) is zero, and all the spin remain down 
(region A in Fig. 3.2). 

For h > 2J — A, and A < J, pi is 1 whenever po is nonzero. Then from Eq. A.l, 
P*{h) becomes 1. Thus, for A < J, we get a rectangular loop and the system changes 
from all spins down to all spins up state in a single big avalanche. 

For A > J,pi—po equals J/ A and is independent of h, in the range 2J — A < h < A. 
Thus P*{h) is a linear function of h in this range, increasing from to 1. 

Defining 

If h 2J\ 

we obtain the expression for P* as 

r for e < 0, 

= \ TTT/A for < e < 1 - J /A, (A.3) 

i 1 for e > 1 - J/ A. 

Using Eq. (3.6), the expression for Qo is, 

Qo = {l-pi)-{p2-pi)P''{h). (A.4) 

The generating function Q{x) obtained from the self-consistent equation [Eq. (3.9)] is, 

= + (A.5) 

1 - x{pi - Po) 
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and the generating function G{x\h) given by Eq.( 3.12) becomes , 

Gix\h) = X {[Q{x)]^^{2J -h) + 2P*[Q{x)]^{-h) + P*2^(-2J - h)} . (A.6) 

Now if A > 2 J, and - A + 2 J < /i < A - 2 J (region B in Fig. 3.2), 

{P2 -Pi) = {pi -Po) = ^/A, 

(j){2J - mJ - h) = ^ for m = 0,l,2; 

and P* + g„ = 1 - J/A. 

1 - (J/A) 



Thus 



Q{x) 



1 - (J/A) 



p. 



and 



(A.7) 



(A.8) 



Expanding G{x\h) in powers of x, we get the probability distribution of avalanches in 
region B given by Eq. (3.14) of sec. 3.2.1. 

In the region C, p2 saturates to value 1, (f>{—2J—h) becomes zero and {p2—pi) becomes 
(1 — J/A — e). Thus we get, 

^ (1 - J/ A - ef 



In terms of P* and Qo we get 

Q{x) 

and 



(1-J/A) • 

)o + xP*[l - 2( J/A) - e] 
1 - i J/A)x 



G{x\h) = ^{[p^+Q{x)r-p^'}. 



Expanding G{x\h) in powers of x we get , in region C 

Giih) = ^[ip^+Qor-p 

and 



J 



Gs{h) = [A[s + A',][^-j , for s>2. 

Here A2 and B2 have no dependence on s but are explicit functions of h 

1 



(A.9) 

(A.IO) 
(A.ll) 

(A. 12) 
(A.13) 



2A 



(J/A 



1 / 2J h\ 
^4(J/A)3 A:~ aJ 



*2 



A 



*2 
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Integrating over h from — oo to oo we get the integrated avalanche distribution Dg, 



G,{h)dh, 



(A.14) 



where 



1 



(1 - J/A)2 
and, for s > 2, 



with 



A2 
B2 



(A.15) 
(A.16) 



1 



A.2 Avalanche distribution on a three coordinated Bethe 
lattice 

For z = 3, the self-consistent equation, for P*{h) [Eq. (2.9)] is a quadratic equation, 

[{P2 - Pi) - {Pi - Po)]P*{hY + [2(pi - po) - l]P*{h) +po = 0. (A.17) 

For the rectangular distribution, the coefficient of P*^ is zero for a range of /^-values, and 
P*{h) is still a piece wise linear function of h 

r for e < 0, 

P"ih) = for < e < 1 - 2( J/A), (A.18) 
[l fore > 1 - 2(J/A), 

where e is defined as. 

The self-consistent equation for Q{x) [Eq. (3.9)] becomes, 

4pi - Po) [Q{x)f + [2xP*{p2 -Pi) -I] Q{x) + xP*\ps - + go = 0, (A.20) 
where Qo is obtained [Eq. (3.6)] as 

Qo = {l-Pi)- 2{p2 - Pi)P* + [{P2 -Pi)-{P3- P2)] P*', (A.21) 
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and the expression for G{x\h) [Eq. (3.12)] becomes, 



G{x\h) = X [[Q{x)f(j){3J -h) + 3[Q{x)]^P*(j){J - h) 

+3[g(x)]P*V(-'/ -^) + P*^</'(-3J - h)] . (A.22) 



Now in the region B, 

(P3 - P2) = {P2 - Pi) = {Pl - Po) = ^/A, 

(f>{3J — 2mJ — h) = -^-^ for m = to 3, 
and P* + g„ = 1 - J/ A. 

Solving Eq. (A. 20) and choosing the root which is well behaved for x near 0, we get 



1 - Jl-4{J/A)x{P* + Qo) 

and the expression for G{x\h) [Eq. (A.22)] becomes 

Gix\h) = ^[P* + Qix)f . (A.24) 

Expanding G{x) in power series of x, we obtain the Eq. (3.19) of sec. 3.2.2. 

In the region C, p3 saturates to the value 1, </>(— 3 J — h) becomes zero and [p^ — P2) is 
no longer independent of h. Substituting the appropriate expressions, we find that 



p. . 1 - - 4(^/A)x[(l - 3(J/A) - .) + (P- + Qo)] 

and 

We note that the term inside the radical sign in Q{x), and also in G{x\h), is a simple 
linear function of x. It is thus straightforward to expand it in powers of x using binomial 
expansion. This gives us the Eq. (3.23) of sec. 3.2.2. 
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